Modulational instability, quantum breathers and two-breathers in a frustrated ferromagnetic spin lattice under an external magnetic field*

Project supported by the National Natural Science Foundation of China (Grant No. 11604121), the Scientific Research Fund of Hunan Provincial Education Department, China (Grant Nos. 16B210 and 16A170), and the Natural Science Fund Project of Jishou University, China (Grant No. jdx17036).

Su Wanhan, Xie Jiayu, Wu Tianle, Tang Bing
College of Physics, Mechanical and Electrical Engineering, Jishou University, Jishou 416000, China

 

† Corresponding author. E-mail: bingtangphy@jsu.edu.cn

Project supported by the National Natural Science Foundation of China (Grant No. 11604121), the Scientific Research Fund of Hunan Provincial Education Department, China (Grant Nos. 16B210 and 16A170), and the Natural Science Fund Project of Jishou University, China (Grant No. jdx17036).

Abstract

The modulational instability, quantum breathers and two-breathers in a frustrated easy-axis ferromagnetic zig-zag chain under an external magnetic field are investigated within the Hartree approximation. By means of a linear stability analysis, we analytically study the discrete modulational instability and analyze the effect of the frustration strength on the discrete modulational instability region. Using the results from the discrete modulational instability analysis, the presence conditions of those stationary bright type localized solutions are presented. On the other hand, we obtain the analytical expressions for the stationary bright localized solutions and analyze the effect of the frustration on their emergence conditions. By taking advantage of these bright type single-magnon bound wave functions obtained, quantum breather states in the present frustrated ferromagnetic zig-zag lattice are constructed. What is more, the analytical forms for quantum two-breather states are also obtained. In particular, the energy level formulas of quantum breathers and two-breathers are derived.

1. Introduction

In recent years, many researchers have paid close attention to studies on solid based quantum computing. In order to realize the transfer and storage of a quantum-state, the spin chain system may be employed in quantum data bus as well as quantum memory.[15] Physically, Heisenberg spin chains are discrete lattice systems, which contain the intrinsic nonlinearity caused by the spin–spin exchange coupling, on-site anisotropy, Dzyaloshinskii–Moriya interaction, octupole–dipole coupling, and so on.[6] It is well known that the nonlinearity lattice system can support a quite interesting nonlinear excitation, i.e., discrete breather (also called intrinsic localized mode).[7] Such excitations are time-periodic and space-localized. The energy of the discrete breather is exponentially localized in space for the short-range interaction between particles. Experimentally, discrete breathers have been identified in antiferromagnetic materials.[8] In this sense, it is quite natural to consider the application of discrete breathers in quantum-information storage.

Until now, a lot of works on (spin) discrete breathers in classical ferromagnetic or antiferromagnetic spin chains have been reported.[922] Theoretically, the spin chain is a pure quantum system without any classical analogue. Classical treatment may be reasonable only when the spin value is very large. It is undeniable that those classical works discovered some important properties of magnetic discrete breather modes. For example, an isotropic spin chain cannot produce the discrete breather.[19] Once the on-site or exchange interaction anisotropy is introduced, discrete breather modes are possible in spin chain systems.[23] Furthermore, Lai et al. demonstrated that discrete breather mode can exist in the classical isotropic ferromagnetic chain involving the second neighbor coupling.[10] Their results seem to show that the anisotropy[24] or long-range coupling[25] plays a vital role in the existence of discrete breathers in isotropic spin chains. Since the spin chain is a quantum system, it is meaningful to investigate intrinsic localized spin wave modes in spin chains within the framework of quantum mechanics. The quantum counterpart of the intrinsic localized mode is termed the quantum breather.[26]

In order to study quantum breathers in spin chains, different approaches have been employed in the relevant published articles.[2733] The first step in all studies is to quantize the spin chain model by using the Holstein–Primakoff or the Dyson–Maleev transformation. Thus, the corresponding quantized spin chain system can be described via an expanded Bose–Hubbard (EBH) model with a certain number of bosons. When the number of bosons is small, one can consider a numerical diagonalization approach combined with perturbation theory. With this approach, Djoufack and co-workers calculated numerically the energy spectra of some ferromagnetic chains.[27,28] Their results showed that these ferromagnetic chains can support 2, 4, and 6-boson quantum breather states. On the other hand, we have developed an analytical approach to study the spin chain system including a large number of bosons. The basic thought for this analytical approach is the application of the (time-dependent) Hartree approximation. By applying this approach, we have constructed the analytical forms for quantum breather states in some ferromagnetic spin lattice chains.[2932] Very recently, Djoufack et al.[33] studied quantum breathers in a weak Heisenberg ferromagnetic spin lattice by means of our approach. Some new features for quantum breathers were found in this weak ferromagnetic spin lattice. They also checked the reliability of the analytical results by making use of a split-step Fourier numerical method.

In this study, we will consider a frustrated zig-zag ferromagnetic spin lattice. In recent years, strongly frustrated low-dimensional magnets have attracted much attention.[34] A quite interesting class of frustrated spin systems can be viewed as chain compounds, which are constitutive of edge-sharing CuO4 components.[3537] In recent years, all kinds of those copper oxides have been composited and shown to reveal particular physical features.[3840] The aim of the present research is to construct quantum breather and two-breather states by making use of our analytical scheme. So far, little attention has been paid to quantum two-breather state in spin lattice system. Hence, the study of quantum two-breather in the frustrated ferromagnetic chain is a new attempt for us. In particular, a key step in our research scheme is that quantum breather states are constructed by using the stationary localized single-quanta wave functions. In theory, the emergence of the bright type stationary localized solution is compactly linked to the discrete modulational instability of spin plane waves.[41] Based on this consideration, the discrete modulational instability of the spin plane wave with the invariable amplitude will be analyzed.

The framework of this article is as follows. In Section 2, we present a zig-zag ferromagnetic lattice under a external magnetic field and quantize the model Hamiltonian. In Section 3, we obtain the equation of motion corresponding to the single-magnon wave functions within the Hartree approximation. In Section 4, the discrete modulational instability of the spin plane wave with invariable amplitude is investigated via employing a linear stability analysis. In Section 5, analytical expressions of quantum breather states are obtained. In Section 6, we obtain the quantum two-breather state in the present magnetic lattice. The last section is devoted to summarizing our work.

2. The frustrated ferromagnetic lattice model and its quantization

In this research, a frustrated zig-zag spin lattice with anisotropic ferromagnetic first neighbor and antiferromagnetic second neighbor exchanges is considered. When the frustrated zig-zag lattice is in an external magnetic field, the corresponding Hamiltonian reads

where corresponds to the μ component for the local spin operator on the j-th site, J1 (< 0) is the first neighbor coupling parameter, J2 (> 0) is the second neighbor coupling parameter, Δ1 is the anisotropy constant of the first neighbor exchange interaction, Δ2 is the anisotropy constant of the second neighbor exchange interaction, and the additional magnetic field B is supposed to be applied along the ez direction.

For the sake of convenience, we shall restrict our research to the system only with the anisotropy of the second-neighbor interaction (i.e., Δ2 = 1). In this case, Hamiltonian (1) takes the following form:

with Δ1 = Δ > 1. In Hamiltonian (2), we have set J1 = −J and J2 = λJ. Here, λ ≡ |J2/J1| is known as the frustration parameter.[42] Dmitriev and Krivnov[43] have studied magnetic solitons in the classical frustrated ferromagnetic zig-zag spin lattice. By that study, they have analyzed the behavior of soliton in the shadow of the ground state phase transition from the ferromagnetic to the spiral phase. In our work, we assume that the zig-zag spin chain system is always in the pure ferromagnetic phase. When the additional magnetic field is sufficiently strong, this assumption should be reasonable.[44]

Taking advantage of the raising and lowering operators , then the Hamiltonian of the zig-zag spin chain becomes

Here, the operators , , and obey the following two commutation relations: and . To quantize the frustrated zig-zag spin lattice that we consider, we may make use of the Dyson–Maleev transformation,[45,46] i.e.,
where the Planck constant ħ has been fixed as the unit, S represents the value of spin length, and and cj correspond to Bose creation and annihilation operators, respectively. The substitution of Eqs. (4)–(6) into Eq. (2) leads to a bosonized Hamiltonian that reads
where ω0 = BB + 2(Δλ)J. In the above bosonized Hamiltonian, the zero-point energy of the present frustrated zig-zag ferromagnetic chain has been set to be zero. Physically, this is equivalent to the zero point for the energy of the ferromagnetic spin chain being reset.

3. Basic physical equations

By adopting the Dyson–Maleev transformation, the frustrated zig-zag ferromagnetic lattice has converted into a system of (local) magnon gas in a zig-zag lattice, which is described via an EBH model with the second neighbor coupling. Here, the Schrödinger picture shall be utilized to analyze the quantum dynamics in this EBH model. Physically, the temporal evolution law for the boson gas system state vector |Ψ(t)⟩ has to obey the (linear) Schrödinger equation

It is obvious that the expanded Bose–Hubbard Hamiltonian (7) definitely commutes with magnon-number operator (the corresponding eigenvalue is n). It is easy to conclude that the overall number of magnons has to be conserved. Hence, we shall adopt a Fock representation to depict the quantum state of the present lattice system. In the Fock representation, a general n-magnon state vector of the system may be written as[47]
where |0⟩ = |0⟩1 |0⟩2 ⋯ |0⟩f corresponds to a quantum vacuum state. Physically, the term Fn(j1,j2, . . ., jn,t) is understood as a n-magnon wave function satisfying the normalization condition[48]
Inserting the expanded Bose–Hubbard Hamiltonian (7) and the quantum state vector (9) into Eq. (8) and utilizing the Bose commutation relation , we obtain the equation of motion of the n magnon wave function Fn (j1, j2, . . ., jn,t), namely,
Actually, equation (11) is a time-dependent Schrödinger equation corresponding to one system of n-magnon in the zig-zag lattice of f sites. Furthermore, we notice that the interaction of magnonmagnon is the Kronecker δ-function. Technically, it is quite complicated to solve the Schrödinger equation (11) exactly for the case of large n and f ≥ 3. In this study, the time-dependent Hartree approximation is to be considered. The approximation approach is feasible when a given number of magnons is sufficiently large.

Under the Hartree approximation, those n-magnon wave functions can take a product of the form[49]

Here, Φn,jk (t) is a single-magnon wave function, where jk = 1,2, . . ., f, and k = 1,2,. . .,n labels all magnons. Taking into consideration that such single-magnon wave functions do not depend on the index k, so they may be directly described via Φn,j(t), where j = 1,2, . . ., f.

By the use of the n-magnon Hartree wave functions in Eq. (12), the n-magnon state vector (8) can be expressed as[50]

Then, the normalization condition becomes
By taking advantage of the variational principle, one can obtain equations of motion for those single-magnon wave functions Φj (j = 1,2,. . ., n)
In the above equation, we have abandoned the suffix n and substituted Φj for Φn,j.

4. Discrete modulational instability

Theoretically, the practicability of bright type localized solutions can be deduced by making use of the modulational instability analysis. Here, we intend to study the discrete modulational instability for Eq. (15) through a linear stability analysis.[51] We begin to analyze the modulational instability via seeking the discrete plane spin wave solution to Eq. (15), which reads

where ϕ0 corresponds to a constant amplitude, q stands for the wave number of the discrete plane spin wave, a is the system lattice constant, and ω denotes the radian frequency of the discrete plane spin wave. The insertion of the discrete plane wave solution (16) into Eq. (15) leads to a nonlinear dispersion relationship
Next, one can add a small perturbation, δΦj, in Eq. (16), and suppose that the solution of Eq. (15) is
Substituting the perturbation solution (18) in Eq. (15) and abandoning those nonlinear terms on the perturbations δϕj and , one can obtain a linear differential equation
In light of the process of the linear stability analysis, we should suppose that equation (19) has a general solution,[52] namely,
where Q and Ω correspond to the wave number and the radian frequency of the small perturbation fluctuation, respectively. What is more, both χ and ς are arbitrary real variables. By inserting the general solution (20) into Eq. (19), one can have
The emergence circumstance for Eq. (21) with the untrivial solution is given via the second order equation for the perturbation circular frequency Ω, which is
with
By doing a simple mathematical derivation, equation (22) becomes
While the right-side term in Eq. (24) is negative, it has two plural solutions for the perturbation frequency Ω. In this case, it is evident that the discrete plane wave (18) is modulationally unstable. In general, the growth rate (gain) of the modulational instability is tightly linked to the imaginary part of this perturbation frequency Ω.[5356] In our physical scenario, the corresponding instability growth rate can be expressed as
where Im signifies the imaginary part corresponding to the perturbation circular frequency Ω. Only when the discrete modulational instability of the plane spin wave happens, can the bright localized solution to Eq. (15) exist.

For the nonlinear equations, the modulational instability provides a significant mechanism to understand emergence of the bright type localized structures. In Fig. 1, we present the modulational instability region on plane (q,Q) corresponding to the present zig-zag spin lattice with various values of the frustration parameter λ. From Fig. 1, one can find that configuration of the modulational instability region apparently changes with the increase of λ. In particular, a new stability region appears near the Brillouin zone center when λ is large enough, as shown in Figs. 1(c) and 1(d). In Fig. 2, we plot the growth rates of the modulation amplitude for the q = 0 plane spin wave. It is found that the growth (gain) of the instability can be influenced by varying λ. We can see clearly that the q = 0 plan spin wave is not modulationally stable for small Q when λ ≤ 0.24. Based on this fact, it is easy to predict that the stationary bright type localized solution can occur at the Brillouin zone center only when λ is less than a specific value.

Fig. 1. (color online) Modulational instability regions (bright zones) on the plane (q,Q) for the frustrated zig-zag spin chain with different frustration parameters: (a) λ = 0, (b) λ = 0.2, (c) λ = 0.3, and (d) λ = 0.4. The other parameters are chosen as J = 1, S = 20, Δ = 2, n = 30, a = 1, and ϕ0 = 0.1.
Fig. 2. (color online) The growth rate (gain) of the instability versus the perturbation wave number for various values of λ. The main parameters are the same as those in Fig. 1.
5. The stationary bright localized solutions and quantum breather states

With the aim of constructing the quantum breather state, we need to seek the bright stationary localized solution to Eq. (15). Here, we adopt a semidiscrete multiple-scale method[57] to solve Eq. (15) approximately. On the basis of this method, one first makes a transformation Φjεϕj. The next step is to consider a trial solution with the modulated wave of the form

where t1 = εt and t2 = ε2t correspond to two multiple-scale time variables, and x1 = ε ja is a new spatial variable; θj = qjaωt signifies the phase of carrier (spin) waves. The key thought for this semidiscrete approximation is a hypothesis that the envelope function ϕ is considered to be continuous meanwhile the phase θ is considered as a discrete function. After executing a simple calculation, we have
It is apparent that equation (27) is the magnon dispersion relationship in the present frustrated zig-zag lattice. In Fig. 3, we display the magnon dispersion curves for different values of the frustration parameter. From Fig. 3, we can see that the configuration of the magnon dispersion curve significantly varies as λ increases. When λ is less than 0.25, the magnon frequency exists a minimum point at the first Brillouin-zone center q = 0 meanwhile a maximum value point occurs at the Brillouin-zone boundary q = π/a. However, q = 0 becomes a local maximum point once λ is greater than 0.25. What is more, a new maximum then appears at wave number q = q0 = cos−1[1/(4λ)].

Fig. 3. (color online) The magnon dispersion curves in the first Brillouin zone for the present zig-zag ferromagnetic chain with different values of the frustration parameter. The corresponding parameters are set to J = 1, a = 1, g = 1, μB = 1, S = 20, B = 50, and Δ = 2.

Applying the dispersion relationship (27), one can calculate the group velocity for carrier waves

Actually, equation (28) contains the relation between the spatial and temporal derivatives for ψ. Based on this fact, we infer that the function ψ should take on travelling waves with the following form:
where z1 = x1Vgt1 is also a spatial scale. Taking advantages of Eqs. (28) and (30), and introducing a transformation ψ = U/ε and z = z1/ε, one can recast Eq. (28) into the form
with
We notice that equation (32) is exactly the famed (cubic) nonlinear Schrödinger equation (NLSE).

Considering that our task is to find the bound state solution corresponding to the single-magnon wave function in the present physical context, hence we will only consider the case of sign(PQ) > 0. For sign(PQ) > 0, the bright soliton solution to NLSE (41) has been constructed by some proven techniques, such as the inverse scattering method,[58] the Darboux transformation,[5963] and the Hirota bilinear method.[6467] Generally, the bright soliton solution to the cubic NLSE takes the following form:

Here, A0 denotes the bright soliton amplitude and z0 is the corresponding initial center position.

Since our idea is to use the bright stationary localized solution to construct quantum breather states, we shall focus on the case of Vg = 0. It should be mentioned that such wave numbers corresponding to Vg = 0 appear at the extreme points of the magnon frequency ω(q). It has been found that extreme points occur at different wave numbers for λ < 0.25 and λ > 0.25. Therefore, we will discuss these two different cases respectively.

First, we consider the case of the frustration parameter λ < 0.25. In this case, the extreme points of the magnon frequency are q = 0 and q = π/a. However, one can find that sign(PQ) is negative at the boundary of the Brillouin zone q = π/a. Hence, the wave number q = π/a has to be discarded. Thus, we need to focus our attention on the Brillouin zone center q = 0. At wave number q = 0, one can obtain

With the results obtained in the above text, we can obtain the stationary localized single-magnon wave function on lattice site j, namely,
with C1 = (Q/2P)q=0 and Ω1 = ωmin − (Δ − 1)(n − 1)Ja2C1/4. Here, Ω1 is the eigenfrequency of the localized solution, which is below the base of the magnon frequency band. Using these stationary localized single-magnon wave functions, it is easy to construct the following Hartree product eigenstate:
It is obvious that the above approximate Hartree n magnon state (40) is a localized quantum (stationary) state. Utilizing Eq. (40), one can calculate the average number of magnons on lattice site j, which is
From Eq. (41), we conclude that all magnons are almost localized on those lattice sites near to the central position j = j0 for the quantum state (40). Theoretically, when the quantum breather state is present at a lattice system, the corresponding result is a quantum localized (in space) excitation, which is completely analogous to the classical counterpart of the intrinsic localized mode.[67] Therefore, the localized Hartree n magnon state (40) is the quantum breather state. In Fig. 4, we exhibit the temporal evolution for the average number of magnons on the j-th site corresponding to the quantum breather state (40). From Fig. 4, we can see that most magnons are localized on those sites nearby site j = j0 and the spatial distribution of magnons still retains over time.

Fig. 4. (color online) Temporal evolution of the average number of magnons on the j-th site ⟨nj(t)⟩(H) corresponding to the quantum breather state (40). The relevant parameters are fixed as J = 1, a = 1, S = 30, Δ = 1.15, n = 20, λ = 0.2, and j0 = 300.

Furthermore, it is easy to derive the Hartree energy formula of the quantum breather state (40), which is

The above equation is evidently a discrete energy level formula of the present zig-zag ferromagnetic system. This manifests that the energy of the zig-zag lattice system for the quantum breather state (40) is quantized.

Now, we turn to the case of λ > 0.25. In this case, the extreme points of the magnon frequency are q = 0, q = q0 = π/a − cos−1(1/4γ), and q = π/a. As was pointed out above, the wave number q = π/a has to be abandoned due to sign(PQ) < 0. Furthermore, sign(PQ) is also negative at the center of the Brillouin zone. So, we only need to consider the wave number q0. At q = q0, we have

When λ satisfies the relationship , sign(PQ) must be positive. Then, we also obtain a stationary localized single-magnon wave function, which has the following form:
with C2 = (Q/2P)q=q0 and
Ω2 is the corresponding eigenfrequency that is located below the base of the magnon frequency band. Similarly, with the localized solution (46), we can also construct the following quantum breather state:

For the quantum breather state (47), we can easily calculate the corresponding Hartree energy, which is

6. Quantum two-breather states

In this section, we shall use the same method to construct the quantum two-breather state in the present ferromagnetic zig-zag lattice system. Since two groups of magnons with overlap are fully similar to a group of magnons in the Hartree approximation, we only need to consider the two independent groups of magnons. Despite magnon–magnon interaction in the same group, each magnon may be assumed to have the single-magnon wave function with the same form. Moreover, magnons in distinct groups act differently and thus, their wave functions have different forms. Based on these considerations, we may construct a quantum two-breather state for the ferromagnetic zig-zag lattice system containing n = n1 + n2 magnons, where n1 and n2 magnons are bound together, respectively. In the Hartree approximation, the total wave function for the present ferromagnetic zig-zag lattice system can be written as

In Eq. (49), the summation extends all over the set {N}, which implies all possible arrangements of [1,2,. . .,n1 + n2] with the grouping of magnons into [1,2,. . .,n1] and [n1 + 1,n1 + 2,. . .,n1 + n2] unchanged. Since the many-magnon wave function Fn1,n2 is symmetric for the index jk, the summation is present. By utilizing the variational approach, one can obtain
In the above calculation, we have considered that those single-magnon wave functions Φn1,j and Φn2,j have to be well-splitted. From the results obtained in the above section, the solution to Eq. (50) depends on the value of λ. Here, we first use the case of λ < 1/4 as an example to construct the quantum two-breather state in the present zig-zag ferromagnetic lattice. For the case of λ < 1/4, the bright stationary localized solutions to Eq. (50) are given by
with
Using these bright stationary localized solutions, it is easy to construct an n magnon Hartree state, which reads
Taking into consideration the normalization condition(H)Ψn(t)|Ψn(t)⟩(H) = 1, we can determine the normalization constant Nn1,n2, that is
With Eqs. (51) and (55), the n magnon Hartree state (54) can be written as
Physically, the classical two-discrete breather is the bound state of two single discrete breathers.[69] When the quantum two-breather state superposes, the result is a localized excitation that is fully analogous to the classical two-discrete breather.[70] Hence, it is not difficult to conclude that the n magnon Hartree state (56) is a quantum spin two-breather state. Furthermore, in Fig. 5, we show the temporal evolution of the average number of magnons on the j-th site corresponding to the quantum two-breather state (56). From Fig. 5, one can observe that two independent groups of magnons are localized on some sites in two different areas.

Fig. 5. (color online) Temporal evolution of the average number of magnons on the j-th site ⟨nj(t)⟩(H) corresponding to the quantum two-breather state (49). Those relevant parameters are fixed as J = 1, a = 1, S = 40, Δ = 1.15, n1 = 20, n2 = 30, λ = 0.2, , and .

Furthermore, it is easy to calculate the Hartree energy of the quantum two-breather state (56). The result is

which clearly indicates that the energy of the present ferromagnetic zig-zag lattice system is quantized when the quantum two-breather state (56) emerges.

For the case of λ > 1/4, we can also construct the quantum two-breather state in the same way. The corresponding quantum two-breather state has the following form:

where
When the system is in the quantum two-breather state (58), the corresponding Hartree energy is

7. Conclusion

In this research, we have constructed quantum spin breathers and two-breathers in an easy-axis ferromagnetic zig-zag chain in the Hartree approximation. Because the bright stationary localized single-magnon wave functions are applied to constructing quantum spin breather state within the Hartree approximation, the (discrete) modulational instability for the corresponding equation of motion was analyzed. We found that the configuration of the (discrete) modulational instability area apparently varies with the increase of λ. According to the result from the modulational instability analysis, we forecasted that the stationary bright type localized solution can emerge at the Brillouin zone center only when the frustration parameter is smaller than a specific value. To seek the analytical expression for the stationary bright type localized solution, a semidiscrete multiple-scale approach was employed to solve the equation of motion analytically. When λ is smaller than 0.25, a stationary localized bright localized solution appears at the Brillouin zone center. However, if λ > 0.25, then the stationary bright type localized solution occurs at wave number q = q0. Using these stationary bound single-magnon wave functions, we first constructed quantum spin breather states for two different cases, i.e., λ < 0.25 and λ > 0.25. For these quantum breather states obtained, we calculated the corresponding Hartree energy of the present zig-zag ferromagnetic, which suggests the energy of our quantum spin breathers is quantized. Furthermore, we constructed analytical forms of quantum spin two-breather states in the present frustrated ferromagnetic zig-zag in the Hartree approximation. With progress of the experimental technology,[7178] our theoretical predictions are expected to be identified in experiments.

Reference
[1] Bose S 2003 Phys. Rev. Lett. 91 207901
[2] Zhang G F Li S S 2005 Phys. Rev. 72 034302
[3] Shi T Li Y Song Z Sun C P 2005 Phys. Rev. 71 032309
[4] Yung M H Benjamin S C Bose S 2006 Phys. Rev. Lett. 96 220501
[5] Lu J Zhou L Kuang L M Sun C P 2009 Phys. Rev. 79 016606
[6] Zakeri K 2014 Phys. Rep. 545 47
[7] Flach S Gorbach A V 2008 Phys. Rep. 467 1
[8] Sato M Sievers A J 2004 Nature 432 486
[9] Lai R Kiselev S A Sievers A J 1996 Phys. Rev. 54 R12665
[10] Lai R Kiselev S A Sievers A J 1997 Phys. Rev. 56 5345
[11] Lai R Sievers A J 1998 Phys. Rev. 57 3433
[12] Kim S W Kim S 2002 Phys. Rev. 66 212408
[13] Rakhmanova V Shchegrov A V 1998 Phys. Rev. 57 R14012
[14] Lai R Sievers A J 1998 Phys. Rev. Lett. 81 1937
[15] Huang G Zhang S Hu B 1998 Phys. Rev. 58 9194
[16] Huang G Velarde M G Zhu S 1997 Phys. Rev. 55 336
[17] Huang G Xu Z Xu W 1993 J. Phys. Soc. Jpn. 62 3231
[18] Sato M Sievers A J 2005 Phys. Rev. 71 214306
[19] Speight J M Sutcliffe P M 2001 J. Phys. A: Math. Gen. 34 10839
[20] Lakshmanan M Subash B Saxena A 2014 Phys. Lett. 378 1119
[21] Khalack J M Zolotaryuk Y Christiansen P L 2003 Chaos 13 683
[22] Fleurov V Zolotaryuk Y Flach S 2001 Phys. Rev. 63 214422
[23] Wallis R F Mills D L Boardman A D 1995 Phys. Rev. 52 R3828
[24] Bai Y H Zhou W P Hou X J Yun G H Bai N 2011 Acta Phys. Sin. 60 056805 in Chinese
[25] Soltani M R Mahdavifar S Mahmoudi M 2016 Chin. Phys. 25 087501
[26] Fleurov V 2003 Chaos 13 676
[27] Djoufack Z I Kenfack-Jiotsa A Nguenang J P Domngang S 2010 J. Phys.: Condens. Matter 22 205502
[28] Djoufack Z I Kenfack-Jiotsa A Nguenang J P 2012 Eur. Phys. J. 85 96
[29] Tang B 2016 Int. J. Theor. Phys. 55 3657
[30] Tang B Li D J Tang Y 2014 Chaos 24 023113
[31] Tang B Li D J Tang Y 2014 Phys. Status Solidi 251 1063
[32] Tang B 2017 Int. J. Theor. Phys. 56 2310
[33] Djoufack Z I Tala-Tebue E Nguenang J P Kenfack-Jiotsa A 2016 Chaos 26 103110
[34] Mikeska H J Kolezhuk A K 2004 Quantum Magnetism Berlin Springer-Verlag
[35] Mizuno Y Tohyama T Maekawa S Osafune T Motoyama N Eisaki H Uchida S 1998 Phys. Rev. 57 5326
[36] Hase M Kuroe H Ozawa K Suzuki O Kitazawa H Kido G Sekine T 2004 Phys. Rev. 70 104426
[37] Malek J Drechsler S L Nitzsche U Rosner H Eschrig H 2008 Phys. Rev. 78 060508(R)
[38] Drechsler S L Volkova O Vasiliev A N Tristan N Richter J Schmitt M Rosner H Málek J Klingeler R Zvyagin A A Büchner B 2007 Phys. Rev. Lett. 98 077202
[39] Capogna L Mayr M Horsch P Raichle M Kremer R K Sofin M Maljuk A Jansen M Keimer B 2005 Phys. Rev. 71 140402(R)
[40] Masuda T Zheludev A Bush A Markina M Vasiliev A 2004 Phys. Rev. Lett. 92 177201
[41] Kivshar Y S Salerno M 1994 Phys. Rev. 49 3543
[42] Dmitriev D V Krivnov V Y 2009 Phys. Rev. 79 054421
[43] Dmitriev D V Krivnov V Y 2010 Phys. Rev. 81 054408
[44] Liu X F Han J R Jiang X F 2010 Acta Phys. Sin. 59 6487 in Chinese
[45] Dyson F J 1956 Phys. Rev. 102 1217
[46] Dyson F J 1956 Phys. Rev. 102 1230
[47] Wright E Eilbeck J C Hays M H Miller P D Scott A C 1993 Physica 69 18
[48] Li D J Mi X W Deng K 2010 Acta Phys. Sin. 59 7344 in Chinese
[49] Lai Y Haus H A 1989 Phys. Rev. 40 844
[50] Wright E M 1991 Phys. Rev. 43 3836
[51] Tchameu D T Tchawoua C T Motcheyo A B 2015 Phys. Lett. 379 2984
[52] Li P Wang L Kong L Q Wang X Xie Z Y 2018 Appl. Math. Lett. 85 110
[53] Cai L Y Wang X Wang L Li M Liu Y Shi Y Y 2017 Nonlinear Dynam. 90 2221
[54] Wang L Zhang L L Zhu Y J Qi F H Wang P Guo R 2016 Commun. Nonlinear Sci. Numer. Simulat. 40 216
[55] Pei S X Xu H Sun T T Li J H 2018 Acta Phys. Sin. 67 054203 in Chinese
[56] Wang L Zhu Y J Wang Z Q Xu T Qi F H Xue Y S 2016 J. Phys. Soc. Jpn. 85 024001
[57] Liu Y Tang Y 2008 Chin. Phys. 17 3841
[58] Ablowitz M Segur H 1985 Solitons and the Inverse Scattering Transform Philadephia SIAM
[59] X Lin F 2016 Commun. Nonlinear Sci. Numer. Simulat. 32 241
[60] Wang L Zhu Y J Wang Z Z Qi F H Guo R 2016 Commun. Nonlinear Sci. Numer. Simulat. 33 218
[61] Wang L Li X Qi F H Zhang L L 2015 Ann. Phys. 359 97
[62] Wang L Zhang J H Liu C Li M Qi F H 2016 Phys. Rev. 93 062217
[63] Wang L Li M Qi F H Xu T 2015 Phys. Plasmas 22 032308
[64] Yang C Li W Yu W Liu M Zhang Y Ma G Lei M Liu W 2018 Nonlinear Dyn. 92 203
[65] Li W Ma G Yu W Zhang Y Liu M Yang C Liu W 2018 Chin. Phys. 27 030504
[66] Liu W Yu W Yang C Liu M Zhang Y Lei M 2017 Nonlinear Dyn. 89 2933
[67] Liu W Yang C Liu M Yu W Zhang Y Lei M 2017 Phys. Rev. 96 042201
[68] Pinto R A Flach S 2008 Phys. Rev. 77 024308
[69] Kivshar Y S Champneys A R Cai D Bishop A R 1998 Phys. Rev. 58 5423
[70] Tang B 2017 Commun. Nonlinear Sci. Numer. Simulat. 48 361
[71] Liu W Liu M Yin J Chen H Lu W Fang S Teng H Lei M Yan P Wei Z 2018 Nanoscale 10 7971
[72] Zhang H Tang D Y Zhao L M Wu X 2009 Phys. Rev. 80 052302
[73] Chen Y Wu M Tang P Chen S Du J Jiang G Li Y Zhao C Zhang H Wen S 2014 Laser Phys. Lett. 11 055101
[74] Ponraj J S Xu Z Q Dhanabalan S C Mu H Wang Y Yuan J Li P Thakur S Ashrafi M Mccoubrey K Zhang Y Li S Zhang H Bao Q 2016 Nanotechnology 18 462001
[75] Guo Z Zhang H Lu S Wang Z Tang S Shao J Sun Z Xie H Wang H Yu X Chu P K 2015 Adv. Funct. Mater. 25 6996
[76] Chen S Miao L Chen X Chen Y Zhao C Datta S Li Y Bao Q Zhang H Liu Y Wen S Fan D 2015 Adv. Opt. Mater. 3 1769
[77] Zhao L M Tang D Y Zhang H Cheng T H Tam H Y Lu C 2007 Opt. Lett. 32 1806
[78] Dai J Zhou P Wang P S Pang F Tim J M Graeme M L Zhang J S Yu W Q 2015 Chin. Phys. 24 127508